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Abstract 

o 

Wc overview the existing OCT work, especially the practical aspects of it. We create a 
^ novel algorithm for 3D OCT segmentation with the goals of speed and/or accuracy while 

^ remaining flexible in the design and implementation for future extensions and improvements. 

' ' The document at this point is a running draft being itcrativcly "developc;d" as a progress 

On report as the work and survey advance. It contains the review and summarization of select 

OCT works, the design and implementation of the OCTMARF experimentation application 
and some results. 
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1 Objectives and Motivation 

Much of the work was done on 2D-segmentation of retinal layers; primarily by processing images 
one slice at a time. Disadvantages of this kind of method include: 

1. relationship between images is not used; 

2. processing time is expensive for OCT-volume data 

As a result, there is more and more 3D segmentation research. Thus, the objectives of this 
work include the review the related work and a proposal of new scheme to segment retinal scans, 
do some PoC experiments. 



2 Literature Review 

This is an overview of the prominent papers in the OCT to comprise the literature review. Some 



of the the specific papers reviewed are Yang et al. (2010); Ghorbel et al. (2011); Fabritius et al. 



(2009); Fuller et al. (2007); Datta et al. (2010a), and others based on their practical aspects. 
What follows are either citations or their summary descriptions in a common algorithmic form 
with the purpose of practical realization in a framework or comparative studies later on. 



2.1 Background Review Papers 



• The original invention of OCT is documented in Huang et al. (1991). 



• Wolfgang Drexler and James G. Fujimoto. State-of-the-art retinal optical coherence to- 
mography. Progress in Retinal and Eye Research, 27:45-88, 2008 



One of the recent student survey papers on the OCT - [Datta et al. (2010a); in 8 IEEE 



pages it discusses the recent (prior and including 2010) approaches to the OCT image 
processing. The authors first discuss in detail the OCT technology from image acquisition 
hardware and software perspectives, applications, visualization, followed by the survey of 
the segmentation literature covering de-noising (speckle noise removal, including media 
filters, wavelets, ZAP, nonlinear, diffusion filters, 7x5 mean filter) and analysis (RFE high 
reflectivity, Markov boundary models, coherence matrix, peak and valley analysis in A- 
scans, adaptive thresholding, various multi-step approaches, active contours, optimal 3D 
graph search, random contour based analysis, support vector machines (SVM), deformable 
fluid-filled models). It is very well illustrated. The authors touch on some of the papers we 



have also briefiy reviewed below, such as Fuller et al. (2007); Fabritius et al. (2009). We 



also have on our list of references intersecting with theirs Huang et al. ( 1991 ); Koozekanani 



et al. (2001); Ishikawa et al. (2005); Fernandez et al. (2005); 



and Fujimoto (2008), which is maybe easier to refer to Datta et al. ( 2010a )'s summary 



Baroni et al. 



(2007b); Drexler 



instead for. The authors also highlight the areas that still need more work in the OCT 
image processing back in 2010, in particular in automation, and some errors in algorithms 
proposed by various authors not being able to cope with either healthy or diseased eyes, 
problems introduced by de-noising (decrease in speed and information loss) Some of their 
concerns have been answered to a degree by newer works published later in 2010 and 



2011 and some of which we reviewed in Yang et al. (2010); Ghorbel et al. (2011). Datta 



et al. (2010a) themselves do not offer newer algorithms in that paper to solve some of the 



highlighted problems. 
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2.2 2D Segmentation Papers 

This is a representative list of the 2D segmentation works that we have not reviewed in much 
detail in here, but provide for reference to the reader. Additional summary reviews may be 
added at a later date. 

• Dara Koozekanani, Kim Boyer, and Cynthia Roberts. Retinal thickness measurements 
from optical coherence tomography using a markov boundary model. IEEE Transactions 
on Medical Imaging, 20(9):900-916, 2001 

• Delia Cabrera Fernandez, Harry M. Salinas, and Carmen A. Puliafito. Automated detec- 
tion of retinal layer structures on optical coherence tomography images. Optics Express, 
13(25):200-216, 2005 

• Hiroshi Ishikawa, Daniel M. Stein, Gadi Wollstein, Siobahn Beaton, James G. Fujimoto, 
and Joel S. Schuman. Macular segmentation with optical coherence tomography. Inves- 
tigative Ophthalmol. Visual Scie., 46:2012-2017, 2005 

• M. Mujat, R. C. Chan, B. Cense, B. H. Park, C. Joo, T. Akkin, T. C. Chen, and J. F. 
de Boer. Retinal nerve fiber layer thickness map determined from optical coherence to- 
mography images. Opt. Express, 13:9480-9491, 2005 

• M. Szkulmowski, M. Wojtkowski, B. Sikorski, T. Bajraszewski, V. J. Srinivasan, A. Szkul- 
mowska, J. J. Kaluzny, James G. Fujimoto, and A. Kowalczyk. Analysis of posterior retinal 
layers in spectral optical coherence tomography images of the normal retina and retinal 
pathologies. J. Biomed. Opt., 12, 2007 

• Maurizio Baroni, Pina Fortunate, and Agostino La Torre. Towards quantitative analysis 
of retinal features in optical coherence tomography. Med. Engin. Phys., 29:432-441, 2007b 

• D. R. Chittajallu, G. Brunner, U. Kurkure, R. P. Yalamanchili, and LA. Kakadiaris. Fuzzy- 
Cuts: A knowledge-driven graph-based method for medical image segmentation. In Pro- 
ceedings of the IEEE Conference on Computer Vision and Pattern Recognition (CVPR), 
2009, pages 715-722. IEEE, August 2009. ISBN 978-1-4244-3991-1. doi: 10.1109/CVPR. 
2009.5206623 

• Delia Cabrera Debuc, Harry M. Salinas, and Sudarshan Ranganathan. Improving im- 
age segmentation performance and quantitative analysis via a computer-aided grading 
methodology for optical coherence tomography retinal image analysis. Journal of Biomed- 
ical Optics, 15(4):484-496, 2010 

• Yazdanpanah Azadeh, Ghassan Hamarneh, and Benjamin R. Smith. Segmentation of 

intra-retinal layers from optical coherence tomography images using an active contour 
approach. IEEE Tran. On Medical Imaging, 30(2):484-496, February 2011. doi: 10.1109/ 
TMI.2010.2087390 

2.3 3D Segmentation Papers 

This section lists majority of the review summaries in the 3D segmentation area. 
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Fuller et al. (2007) (a follow up work on the author's first attempt in Zawadzki et al. ( 2007[ >) 



propose the use of the support vector machines (SVMs) and machine learning to do the 
segmentation of a 3D volume. The approach is interactive and involves comprehensive 
visualization. A 3D volume is constructed and the clinicians interact with it by placing 
control points at the layers of interest and then the system interpolates in between. These 
become them the input to the SVN to learn and further segment retina scans in 3D. In 
this latest approach of the authors a multi-resolution imaging is used of the same area 
of interest in order to combat speckle noise and improve overall precision. The run-times 
are approx. within the range 15-30 minutes depending on the hardware. This is a good 
example of a machine learning segmentation, that claims to do better than neural networks, 
etc. but it also requires interactive user sessions before it can classify reliably the retina 
layers. The summary is in Figure [T| 



1 Load 3D OCT scan and visualize; 

// Interactive step by clinicians 

2 Specify vertices of interest for each layer; 

3 begin 



Location of the voxel in the 3D volume: Pi{xi, yi,Zi); 

Scalar intensity at pji fi] 

PiS neighbors: = {i it 1, j it 1, it 1}; 

Mean intensity value of pj's neighbors: 



Variance around pi: 



'^i ~ |jV| XI [-^^ A 



Gradient around pf. 



10 end 

11 Interpolate; 

12 Perform 4- fold multi-resolution image split; 

13 Interpolate the layer points for each resolution; 

14 Machine-learn the points using kernel-based SVN using features such as gradient, 
intensity, variance, and spacial coordinates; 

15 Classify and measure layer thickness; 



Figure 1: Synthesized from Fuller et al. (2007) 
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• Zongqing Lu, Qingmin Liao, and Fan Yang. A variational approach to automatic segmen- 
tation of RNFL on OCT data sets of the retina. In 16th IEEE International Conference 
on Image Processing, volume 1-6, pages 3309-3312, November 2009 



From roughly the same people Garvin et al. (2008), and Garvin et al. (2009) on 3D OCT 



Review for practical aspects: Fabritius et al. (2009). The authors automatically cover 2 



layers: ILM and RPE, but they do it across the entire 3D volume. It seems they are pretty 
fast in most of their iterative algorithms where the number of iterations is configurable and 
serves as a trade off between quality of the segmentation and speed. Claimed speeds on a 
PC range between 16s to 21s. The ILM classification in particular is worth building upon. 
They use signal intensity variation based segmentation; no denoising on both healthy and 
diseased retinas. The algorithmic summary is in Figure [2] and Figure [3} Their RPE 
and ILM can be done in parallel, so speeding up the overall performance on parallel and 
distributed architectures. 

• V. Kajic et al. Robust segmentation of intraretinal layers in the normal human fovea using 
a novel statistical model based on texture and shape analysis. Optics Express, 18(14), 2010 

• M. Mujat, R. C. Chan, B. Cense, B. H. Park, C. Joo, T. Akkin, T. C. Chen, and J. F. 
de Boer. Retinal nerve fiber layer thickness map determined from optical coherence to- 
mography images. Opt. Express, 13:9480-9491, 2005 



This work of Yang et al. (2010) uses gradient information to get 9 layers; they are 



fast and thorough, but the paper omits a few algorithm details. The first author is from 
Topcon, the maker of the OCT 1000 imaging hardware. Subjects - 38 (19 healthy and 19 
glaucoma), including live and dead tissue. The results are compared to 4 experts' manual 
segmentation. 

The parts of algorithm include two steps - Canny edge detection (customized with 3 
thresholds, graph-based node cost assignment, dynamic programming (during the shortest 
path search). The gradient aspect seems useful and fast. 

The summary in the algorithm-like notation is in Figure [4} 

• Rupsa Datta, S. Aditya, and D. N. Tibrewala. Advancement in OCT and image-processing 
techniques for automated ophthalmic diagnosis. In Proceedings of the 2010 IEEE Students' 
Technology Symposium, pages 26-33. IEEE, April 2010b. ISBN 978-1-4244-5974-2 



The result Ghorbel et al. (2011) of is also on 3D using similar equipment to Yang et al 



(2010) and 8 layers using a different approach. The paper is very detailed in the algorithms 
used. The authors use active contours, Markov random fields, Kalman filters. 

The authors don't seem to claim speed. Denoising is in play. Very detailed literature 
review. They have a lot larger database of 700 images (but only on healthy subjects) and 
test their approach from two separate imaging hardware types: Topcon's 3D OCT 1000 



and Spectralis IIRA+OCT(Heiderlberg) and 5 experts vs. 4 in Yang et al. (2010) 



The methods in Fabritius et al. (2009); Yang et al. (2010) may be more practical 
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// Preprocessing is not done, but can be done if wanted 

// RPE identification 

begin 

// X — width, z — height, y — depth in slices (B-scan) 

Determine coordinates of max. intensity niax{Ix^y{z)) pixels -2max(/)(3^) 2/) across 
volume on (x, y); 

Obtain 2D RPE position matrix 2:max(/)(2^i 2/) = Zrpei{x,y); 
II Takes care of the speckle noise from RNFL 

Using Otsu method, binarize to get GrroriGOUS pixels nicisk Bj-p^ii^x ^ y) for RPE 
position matrix in (j;, y) using top-hat filtering using structuring element 5x5 pixels; 
Replace erroneous pixes with values from nearest neighbors in(a;,y); 
begin 

// Set expected erroneous pixels' value to NaN 
if Brpei{x,y) = then 
I Zrpeiix,y) = NaN; 
end 

Zrpelix , y) =forall the Zrpel {x,y) = NaN do 

Zrpei{x,y) = based on the nearest neighbor value; 
Smooth by moving window median filter 30(x) x 2{y) pixels; 
end 
end 

// 30 initial pixels 

P = {Ix,y{z), z G [Zrpelix, y) - 10, Zrpel{x, y) + 20]}; 
for selected number iterations do 

Extract P pixels around RPE estimation from the original 3D volume {x,y,x); 
Update RPE position based on the maximum intensity; 
Smooth by moving window median filter 40(x) x 2{y) pixels; 
Reduce amount of pixels P around RPE to 20, and then 10 to do 

P = {Ix,y{z),Z £ [Zrpe2{x,y) - 10 , Zrpe2{x, y) + 10]} and 
P = {Ix,y{z),Z £ [Zrpe3ix,y) - 5, Zrpe3{x,y) + 5]}; 

end 

Smooth by moving window median filter 20(x) x 2(y) pixels; 



23 end 



Figure 2: Synthesized RPE Detection from Fabritius et al. (2009) 



3 The Work of Sun et al. at Tsinghua 



3.1 Previous Results 

Sun et al. covered various aspects of the OCT image processing and classification using earth 
mover distance, SVMs, providing for detection of the nacre's layer thickness and other works, e.g. 



Sun and Lei (2009 2010). The most recent work in a 3D algorithm for segmentation (submitted) 
is by Sun and Zhang that explores position determination for ILM, RPE, and IS/OS in a 3D 
volume in an efficient and robust manner by determining intensity change on both sides of the 
boundary simultaneously and then smoothen the recovered surface using 3D intensity difference 



Sun and Zhang (2012). Some segmentation results are shown in Figure 12 and more are within 
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// ILM identification 
1 begin 
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Calculate threshold value and binarize each B-scan slice; 
begin 

foreach B-Scan in data cube do 

Extract first 5 depth (assumed noisy) pixels 
{Noisey{x,z),y e [l,N],xe [l,M],ze [1,5]}; 
Evaluate 5 x M pixels; 

// Bnoisey{x, z) is binarized Noisey{x, z) as each Threshold is 
different 

Threshold = 0.5% of pixels set to after binarization 
{Exli ELi BnoisEyix, z) < 0.005 x 5 x M}; 
end 
end 

Do first ILM estimation via depth position of the first zero value of each binarized 
A-scan; 

for number of iterations do 

Extract 45 pixels around estimated ILM and reprocess to remove erroneous ILM 
pixels; 
begin 

I {Sy{x,z),ye [l,N],x e[l,M],z e [Zum.{x,y) - 15, Zum{x,y) + 
end 

Smooth by moving window medial filters l(x) x 25(y) and 25{y) x l(x); 
Binarize using intensity with the same threshold; 

Re-estimate ILM position based on the first zero depth value of A-scans; 
Repeat with 30 pixels, [Ziim2{x,y) - 3, Ziim2{x,y) + 27], and filter 10(2;) x l{y); 
end 



Figure 3: Synthesized ILM Detection from Fabritius et al. (2009) 



the cited work. 



3.2 Data Description 



Data format description from Sun et al. with the clinical OCT data provided by Shenzhen 
MOPTIM Imaging Technique Co., Ltd. 



• Text files: ASCII-encoded pixel intensity data. These are A-scans, where each line of text 
corresponds to a vertical scan. There are 100 txt files. ( |data.zip^ 



• Binary images: JPEG. 100 jpg files of OCT retinal scans (pic . zip ) representing 100 slides 
of a 3D retina scan at 480x300x100 (96dpi) resolution. 
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Data: 480x512x128 voxels 

Result: 9 segmented boundaries in 16 seconds 
// 2-step segmentation based on gradient: 

1 foreach 4-08x512 A-scan a of the 3D OCT Volume do Obtain gradient information; 

2 begin 



7 
8 
9 
10 



Local customized Canny main edge result with 3 thresholds; 
begin 

I C{i, j) = wi ■ Canny (i, j) + ■ Axial(i, j) + W3 ■ Other s{i, j); 
end 

// To interpolate missing/weak gradient info (over blood vessels, etc) 
Global axial intension gradient; 

Apply the shortest path search to complete and optimize detection; 
begin 

Dynamic programming to optimize the result 

(00 j < 1, j > m 

= < i = n 

y mmm^j-2:j+2 {t(i — + C{i,j) otherwise 

11 end 

12 end 

// Nine boundary detection 

13 foreach B-scan b of the 3D OCT Volume do begin 

14 Preprocess the OCT image with a customized cross-correlation based alignment algorithm to realign the 
A-scans; 

15 Detect the ILM and IS/OS boundaries (the shortest path search maps of the other boundaries can be 
restricted to successively smaller search areas); 

16 Add features edge direction and pixel intensity as the additional cost to the search graph depended on the 
boundary of interest. The edge direction is only considered as black-to-white or white-to-black; 

17 Align the image further the IS/OS boundary; 

18 Detect the OS/RPE and BM/Choroid boundaries; 

19 begin 

20 Accurately detect the OS/RPE and BM/Choroid via a smaller kernel size applied in the Canny edge 
detector to increase the axial detection sensitivity; 

21 The search graph is constructed using only Canny edge and axial intensity gradient strength, and the 
search area is limited below the IS/OS boundary; 

22 end 

23 Detect IPL/INL and NFL/GCL; 

24 Detect GCL/IPL within the NFL/GCL and IPL/INL boundaries; 

25 begin 

26 I Apply the Canny edge detector again within that area to extract the GCL/IPL boundary; 

27 end 

28 Detect the INL/OPL, the dark-to-bright edge between the IPL and the outer nuclear layer; 

29 Detect the ELM similarly via the dark-to-bright edge; 

// Although the detection of a single boundary may utilize the other two pre-detected neighboring 
boundaries to limit the shortest path search area, the intra-retinal boundaries are allowed to 
overlap their neighboring boundaries. 

30 Smooth intra-retinal boundaries (NFL/GCL, GCL/IPL, IPL/INL, ELM, OS/RPE, BM/Choroid) using a 
polynomial curve-fitting based technique; 

31 end 

32 Apply additional smoothing across frames for 3D volumes; 

33 Filter the boundaries for each A-scan location with a ID Gaussian kernel across B-scans; 



Figure 4: High-Level Algorithm for Fast 9-Layer 3D Segmentation by Yang et al. (2010) 



1 Preprocessing; 

2 HRC Detection; 

3 ILM Localization; 

4 Photoreceptor (IS-OS) Segmentation; 

5 Alignment/Clivus Detection; 

6 Inner Layers Segmentation; 



(2011) 



Figure 5: Overall algorithm by Ghorbel et al. 
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1 begin 

2 Crop image to OCT data; 

3 Normalize [0..1]; 

4 Apply non-linear diffusion filter 

5 end 



Weickert et al. 



1 1998 I 



Figure 6: Synthesized Preprocessing from Ghorbel et al. (2011) 



1 begin 

// W is the image width 

2 yo = W/2; 

3 Select pixel {xo,yo) as max response to avg. filter; 

4 Estimate HRC thickness Thrc via profile analysis; 

5 Detect HRC median line; 

6 begin 

7 (1) Vertically smooth image S{x,y) = apply ID Gaussian filter with std. deviation <j = Thrc/'^\ 

8 (2) // Regular median line regardless noise and blood vessels 

9 a = 0.9; 

10 From {xo,yo) column- wrise in both directions deduct median line iteratively from max. output of the 
recursive low-pass filter: max(C(x, i/) = (1 — ce)S{x, y) + aC{x, y it 1)); 

11 (3) Localize HRC contours; 

12 begin 

13 Perform fc-means classification (fc = 3); 

14 From the median line initialize active contour curve X{s) = [x{s), y{s)], s 6 [0, 1]; 

15 begin 
// The active contour is used for regularization purposes 

16 Set a and /3 for contour tension and rigidity; 

17 Compute data fidelity Ej^age{X{s)) as a spatial diffusion of the gradient of an edge map; 

18 Move X within the image's spatial domain to minimize energy E{X): 

E(X) = f i {a\X'(s)\^ + + Ej^,g,{X{s))) ds; 



19 
20 
21 
22 end 







end 



end 



end 



Figure 7: Synthesized HRC Detection from Ghorbel et al. (2011) 



1 begin 



2 


Re-use fc-means classification (fc = 3); 


3 


foreach Image column c scanned top-bottom do begin 






// Pick a boundary pixel p(c, y) if 


4 




switch Pixel do 


5 






case p = first pixel of k = 3 


6 






case p above HRC 


7 






case p = first pixel of k = 2 


8 






Pick p{c, y); 


9 






endsw 


10 






otherwise 


11 






Don't pick p(c, 3/) 


12 






endsw 


13 




endsw 


14 


end 




15 


Select 2 highest gradient pi and pr from the left and right parts of the image; 


16 


Apply edge-tracking algorithm maximizing local gradient starting from pi and pr 


17 


Merge 


resulting curves maximizing the mean gradient; 


18 


Refine ILM boundary by active contour; 


19 


Identify foveola F{xp,yp) where x is maximal on ILM; 


20 end 







(2011) 



Figure 8: Synthesized ILM Localization from Ghorbel et al. 
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1 begin 

2 Apply a peak detector on the image area WRT HRC; 

3 Detect and label maxima in each column to form peak lines; 

4 Init IS/OS junction extraction by selecting the peak line with min distance to HRC; 

5 Iterate over other peak lines WRT the current curve; 

6 begin 

7 Similar mean distance to HRC without overlap; 

8 Apply active contour using detected curve to initialize to get a contiguous regular curve; 

9 end 

10 Detect ONL/IS boundary; 

11 begin 

12 foreach Column n do Apply Kalman filter with a state vector X; 

13 begin 

14 Init from the peak line above IS/OS; 

15 X = [Intensity , Distance); 
// identity matrix 

16 F = I; 

17 Prediction: X{n\n - 1) = FX{n - l\n - 1); 

18 Y(n) = max intensity pixel around X(n\n — 1); otherwise; 
// G(n)adjustedeachiteration 

19 G{n) = weighted Kalman filter gain for error correction; 

20 Update: X{n\n) = X{n\n ~ 1) = G(n){Y{n) - X{n - 1); 

21 Stop: Yin) = for several iterations; 

22 Regularize detective ONL/IS boundary curve by active contour; 

23 if Image is of good quality then 

24 I Apply second Kalman filter to detect OS/RPE boundary similarly to the above; 

25 end 

26 else 

27 Apply fc-means with k = 2 between IS/OS and HRC; 
// Provides first estimation around RPE-ChChap inner boundary 

28 Estimate thickness RPE-ChChap thickness TupE+ChChap', 

29 Perform local analysis of the estimated curve for n; 

30 begin 

31 OS/RPE pixel's y is a first pixel of max gradient below min intensity below the IS/OS 
junction; 

Regularize detected curve by active contour; 

end 

end 



32 
33 
34 
35 
36 
37 end 



end 



end 



Figure 9: Synthesized Photoreceptor (IS-OS) Segmentation from Ghorbel et al. (2011) 



1 begin 



2 
3 
4 

5 end 



Align image by vertical column translation against outside of the RPE+ChCap; 
Define clivus as two highest points xci,yci ^md xcrtVCr on ILM; 
Refine foveola position as max^igtanceiON L/ 1 S, OS/ RPE); 



Figure 10: Synthesized Ahgnment /Chvus Detection from Ghorbel et al. (2011) 
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1 begin 

2 Model noise and spatial pixel interaction using Markov Random Field (MRF); 

3 Classify with Bayesan max a posteriori (MAP) criterion to locate a label config maximizing P based on the 
observed image intensity; 
begin 

// Probability of pixel intensities fs in class i 
/ / cTi — std. deviation 
// fj>i — mean 

1 „ 2^^? . 



6 
7 
8 
9 
10 

11 



12 
13 
14 
15 
16 
17 
18 
19 
20 
21 
22 
23 
24 

25 
26 
27 
28 
29 
30 
31 
32 
33 
34 
35 
36 
37 

38 
39 
40 
41 
42 
43 
44 
45 
46 end 



P(fs\u)s=i) 



/27r<Tj 



// Potts model pixel interactions (f>{ws,wt) in 8-connected sites s and t 
if Ws = i then 

end 

/3 empirically set to balance image data term and regulation term in U ; 
// Energy function 



Uiw\f) = J2 



{fs - ^J,w,)'^ 

2al_ 



+ ln(\/27rCT„^ ) 



Run unsupervised fc-means to init labels and initialize noise params {fj.i,i7i); 
Minimize U using Simulated Annealing repeatedly on the original, but labeled image; 
Determine active contours from the labeled image; 
Detect RNFL; 
begin 

Apply fe-means (fc = 2) to the image top around ILM; 
Select clusters with higher mean intensity; 

Deduce RNFL contours on both sides of foveola with active contours; 
end 

Detect ONL; 
begin 

fc = 3; 
/3 = 2; 

Apply fc-means/MRF segmentation to the foveal region constrained by clivus and ILM — RNFL and 
ONL/IS; 

Deduce the ONL boundary from the labeled image; 
Regularize with active contours; 
end 

Detect INL; 
begin 

fc = 2; 
/3 = 5; 

Apply the above method on the right and left sides from foveola; 
Estimate location of searched boundary as (?/Cii ?/C;/2, ycr/2, J/Cr); 
Force the curve to pass at yp between ILM and OPL/ONL; 
Linear interconnect the 5 points above as the first estimate; 
foreach Column do begin 

Transition pixels between labels 1 and 2 and close to the first estimate arc marked as the new 
boundary; 

Initialize the new boundary estimate as active contour; 

end 
begin 

Traverse column down-up; 

Apply the similar process between the found curve and OPL/ONL; 

end 
end 



end 



Figure 11: Synthesized Inner Layers Segmentation from Ghorbel et al. (2011) 
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Figure 12: Sun et al. Some Segmentation Results 

4 New Algorithm 

This section constructs a new algorithm or a new algorithm framework derived from the reviewed 
literature on OCT and the previous works of authors brought into the OCT domain. 

1. Use gradient, Canny edge detection. 



2. (2D done) Apply Zhang-Suen transform to cleanse the image and skeletonize Zhang and 
|Suen ( |1984 ). Augment the algorithm to 3D covering the nearing planes and gradient 
information. Use binarization before the transform per Fabritius et al. (2009). 



3. Apply 3D discrete wavelet filter (dual tree) to denoise a 3D volume from Selesnick et al. 
( [I503| ). Potentially borrowing from [Kokare et al.| ( pOOSl |2006[ ). 
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4. (stub) GIPSY-based [The GIPSY Research and Development Groupl (|2002-2012[) dis- 



tributed evaluation of OCT images. 



5. Optionally re-use Neural Network (per Baroni et al. (2007a)), 2D CFE filters Haridas 



(2006), and machine learning implementation from MARF (some OCT papers used neu- 



ral networks to classify pixels belonging to different layers.) The MARF Research and 



Development Group (2002-2012) 



6. Optionally re- use Simulated Annealing from Cryptolysis Mokhov et al. (2005-20 12a), as 
e.g. the authors of Ghorbel et al. (2011) use it, as per Figure 
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7. 



Datta et al. (2010a) in their survey pointed out a work of Fernandez (2005) who provided a 
way to detect fluid-filled areas in pathological retinas using deformable models. Separately, 



Miao Song in her master thesis Song (2007) implemented deformable softbody simulation. 



where a 3D sphere of a spring-mass system is subject to inner pressure force opposed to 
by the spring stiffness and other external forces with collision detection. A novel idea is to 
use Song's softbody to fill in the regions of each layer centered at the middles of each layer 
and bounded by the different intensity /gradient value pixels on the boundaries acting as 
collision detection points. Putting enough simulated inner pressure to the softbody object 
would make it fill up the layer accurately enough in real-time to define two layer boundaries 
for most layers for the whole 3D layer in one simulation. 



Datta et al. (2010a) also mentioned an optimal 3D graph search by Garvin et al. in 2008 



- a novel extension of this idea is to use a 3D graph with probabilities and treat the layer 
detection problem as a model-checking problem with probabilities instead that can be 
learned at 3D level with some initial training/learning. For this, there is an open-source 
PRISM - probabihstic model-checking tool iThe PRISM Teaml (|2004-2012|). A model is 



build for each layer with typical voxels of that layer and the corresponding graph; the 
probabilities are assigned based on the local and neighbors averages or medians as well as 
variance and gradient. Once model is build, graph is built and checked against the model. 

9. Export as OBJ file format as well as 3D visualization of the entire segmented volume in 



OpenGL OpenGL Architecture Review Board (1998-2012) to enhance perception of 3D 



layers and "peek" inside the 3D retina. This is a an easy and cost-effective way to visualize 
the results in 3D instead of just each 2D slice of the volume scan. TODO: cite volview 
and other related visualization work 



5 Design and Implementation 



5.1 Language 

Primary language for experiments chosen for now is Java due to its more formal nature, better 
design, memory management. Several convenient frameworks are available for image manip- 
ulation and patter recognition, machine learning in Java as well as easier file management, 
distributed evaluation, web services, and others the author Mokhov is familiar with. Any Java 
program can interface a C++ program and vice versa via the Java Native Interface ( JNI) [Sterns 
(2001-2005); GCJ Contributors (2007). Most critical and interesting things can be converted 



to C++ IFF needed at a later date when the experiments are over. These design decisions and 
implementation are based primarily on Java 6 (1.6.29). 
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// Preprocessing (optional) 
switch preprocessing method do 

case Skeletonize /Thin/ Fill / Contourize 
Optionally binarize per 



Fabritius et al. 



(2009); 



Selesnick et al. 



(2003); 



Use augmented 3D Zhang-Suen transform Zhang and Suen (1984) to cleanse the 
image and skeletonize it; 

5 endsw 

6 case Wavelet 

7 I Use 3D dual-tree wavelet transform filter 

8 endsw 

9 case None 

10 I Skip preprocessing; 

11 endsw 

12 endsw 
// Analysis and Detection 

13 switch detector method do 



14 
15 
16 
17 
18 

19 
20 
21 

22 
23 
24 

25 
26 



Fabritius et al. 



(2009); 



case ILM 
I 3D ILM detection per 
endsw 
case Canny 

Apply the customized Canny edge detection from Yang et al. (2010) on the 
cleansed image; 
endsw 
case Softbody 



Expand softbody Song ( 2007 ) 3D "bubble" bounded by higher intensity pixels 



centered at the layers' middle until the expansion stops; 
Contour of the bubble will defined the layer boundaries; 
endsw 
case PRISM 

// May need to learn the models 
Build a probabilistic graph layer model; 

Treat the problem as a model-checking problem for each pixel in 3D space with 
PRISM 



The PRISM Team 



(2004-2012); 



endsw 

case Simulated Annealing 

endsw 

30 endsw 

31 Visualize with J3D and OpenGL; 

32 ...; 



Figure 13: Building Our Novel Algorithm 
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5.2 Architecture 

• OCTMARFApp - the main application. The tentatively called application stub for experi- 
ments has been created and named OCTMARF. The application's overall structure at the 



moment of creation was based off MARFCATApp (iMokhov' f2010|- , 2011 )) and its predecessor 



WriterldentApp [Mokhov et al.| ( |2009[ ) ; |Mokhov| ( |2008-2012[ ) and others. 



mop 

12j) a: 



The project has been organized into Java packages for various tasks and it incorporated 
an old work on pattern thinning, skeletonization and 2D feature extraction based on 



Zhang-Suen transform Zhang and Suen (1984) that is planned to be extended to 3D 



per the Ideas section. This code was ported, compiles, and runs; it is found under the 
marf . apps . oct . OCTMARF . framework. a2 and its subpackages as of this writing. Alrready 
used in preliminary experiments. 



• See screenshot of Eclipse on the current project's architectural layout overview in Figure 14 



I tP « b ml » I ffi 



ft Packaas E-M]^sr \. Js Heran 
□ F 

□ r 

I I Hi I i, Oct, OCTMARF 

I IJI ■■jrTMARFABB.iava 
. . Ill (S riCTMARFApp 
: coilact-fila5-aiala.pl 
■it.J-* 
II. 



ct, OCTMARF, fram5woi'k,a2,F5atureExtraction 
tt, OCTMARF, fiamenork.aS.FeatureExtiaction.F 
ct, OCTMARF, frameHQrk,a2, Pattern 



□ -ffl P, 
a-ffl Si 



□ IZl [ 

□ ■IZl IDet 



□ a E zer.jav^ 
a-a IPreproce.™.]ava 
El -a Prep,oce55or.)ava 

E ..cct. OCTMARF. f,am* 

ffl ^ ct. OCTMARF. f,am.. 

El-B lProce;;or,]ava 
El B ZhangSuenTiansFgrm,]«va 

i, oct, OCTMARF, OCT 
III I ■;, oct, OCTMARF, OCT, tie 

i.oct.OCTMARF.5torag. 
a a AnyToWAVLoadei.java 
a {il A5canLoacJer.)ava 

□ Be FC.-nper.,ava 
El aF Inm.java 

El a ,cM.taXMLLoader.]ava 

□ ■a PPMDLmper,]ava 

□ a impleLMder,)ava 

□ a e,]ava 
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m 



i riilWIurM.iava ( B ZhargSuenTfari;Hirrri,i JjT] Pattern, lava [^^^^^j^^^j^^.j^^^^^^ ril BirHfEei'iiava ] B PfeEtotesaa'.iavj ") " 



:ayLi3t<Conllgur 



oCon±ig - new Coiifigucatioii [| ; 

□Conf ic 

□Conlig.bSSCountouElng = false 
oConlig.bZSFilllno = false; 



I Coiifigucatioii [I ; 



oConiig - new Conf igucatioii [| ; 
□ Conf ig.bBijiaiisq = false; 
□Conf ig.bSSCountoUEing = false 
□Config.bSSrilllng = true! 




'E: P.Dblem5 f @ Javadoc fg, C.c 



, Searcti i^l Console E3 '-.^ 

<termirated> OCTMARFApp [Java Application] Ciyrograrii FfaUav3\irel,6,Q_IJ3\bin\]avaw,Me<Jan 2, 2Q12 U 
D^iie. [52] 

Filling Pass: 14... Fillingl: Pattern Siae: 4B0S30O 



3 



[19] 



Fillii 



15. 



?. [□] 

Done. [1772] 
Tliiiiniiig Pass: IB 
D^iie. [-157] 
Thinning Pass: 19 
D^ne. [152] 



Thli 



: 20 



Fillingl: Pattern Size: 450 x 300 

Fillingl: Pattern Size: 460 x 300 

. Thinningl: Pattern Size: 480 x 300 

, Thinningl; Pattecii Size: 480 y. 300 

, Thinningl: Pattern Size: 4B0 >: 300 

. Thinningl: Pattern Size: 450 x 300 



Figure 14: Preliminary Project Layout 



5.3 Data Structures 

• Volume ~ a class to contain the 3D OCT volume data. 

• consists 3D arrays: short, double. String representing all the slice of a retina scan 
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— short - primarily for basic loading from the .txt data and manipulation 

— double - for general processing, filtering, etc. as most are floating point numbers 
and many packages work with that, including MARF. 

— String ~ for convenience for debugging and PPM output and printing the values 

5.4 Data Management 

• AScanLoader - a class to load the provided .txt scan data into the internal data structures. 

• PPMDumper - to dump .ppm image format from the internal data structures 

• BMPDumper - to dump .bmp image format from the internal data structures (not imple- 
mented) 



To preliminary tests of loading, basic manipulation, and dumping are successful, as per 
Figure [16] and Figure [16} These simply illustrate that the data are loaded and interpreted 



correctly from the provided data set in Section 3.2 



5B.ppm-19.0 (RGB, 1 layer) 4B0k300 - GIMP 



File Edit Select View Image Layer Colors Tools Filters Windows Help 



Jnjx] 



1^ 



















1} 



[ pn ■'•[ j 100% ■» Background (1,4 MB) 



Figure 15: Test of Loading Provided Data in PPM To Match Provided 



5.5 APIs 
5.5.1 Internal 

The internal API at the moment centers around the detector /analyzer /preprocessor frameworks 
and is invoked by OCTMARFApp. Most of them take the Volume data structure and do something 
with it in either 2D or 3D or both. The API provides for both working on the actual Volume 
instance or its copy if desired to preserve the original. 
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-5B.ppm-18.0 (RGB, 1 layer) 4B0k300 - GIMP 



File Edit Select View Image Layer Colors Tools Filters Windows Help 



Jnjx] 



















1^ 




1 09 , Oj 1 29 , 1 pii -^1 1 100% ■» ^ Click to paint (Ctrl to pick a color) 



Figure 16: Test of Rudimentary Processing via Threshold 

• IDetector - specifies the API ah boundary detectors should implement 

• Detector - generic abstract class that concrete classes may inherit from for convenience 

• FabritiusILMDetector - a specific detector being implemented 

• ZhanSuenTransf orm - an example concrete detector and a processor that relies on the 



implementation of the Zhang and Suen (1984) 



• IPreprocessor - defines API for all preprocessors (e.g. binarizers, denoisers, or filters or 
whatever) to use 

• Preprocessor - a generic class implementing rudimentary API for convenience 

• Binarizer - a concrete simple preprocessor binarizing images with a set threshold 

• IProcessor - API for all kinds of analyzers 

• Configuration - flags and settings for experiment automation 



5.5.2 External 

External API has to do with libraries, frameworks, used or planned to be used for evaluation, 
processing, or connectivity with other platforms and languages. 



• MARF |The MARF Research and Development Group| ( |2002-2012[ ) 

• JAI 
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JNI / JN A |Sterns| ([2001-2005 ) 



• GIPSY The GIPSY Research and Development Group (2002-2012) 



6 Experiments 

Some experiments are in place to test the data processing, the framework's API operation, and 
some actual algorithms. That includes loading of the . txt data and its subsequent dump back as 



an image work (see Figure 15) as well as testing the data structure by rudimentary processing 
of the loaded data's pixels by simply discarding values below certain intensity threshold and 



amplifying the remaining ones to the max as a very basic test (see Figure 16) to ensure the 
given data understanding is correct via Binarizer. 

More interesting preliminary results include global 2D (2D at the moment, per slice, but 
will expand the window to 3D space next) processing that includes in one algorithm some 
subalgorithms that do preprocessing, skeletonizing, or filling/thinning, and active-contouring of 
the images from the Zhang-Suen transform. There are many of the options and configurable 
parameters to try, but already now the outlines of the ILM, IS/OS, RPE can be seen in a very 



draft debug mode in Figure 19 Figure ^Ol Figure 17, Figure \18\ Figure 21 and Figure 22 



Framework accepts Configuration object where multiple configurations with or without 
preprocessing can be tried on the whole volume. 

These can be further refined with parameter tweaking as well as other algorithms (e.g. Canny 
edge detection, graph search, etc.). 

Run-times for the current experiments varied from 29 seconds to 5 minutes depending on 
the configuration for the whole 3D volume. 




Figure 17: 2D Zhang-Suen Skeletonizing near Fovea 
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Figure 18: 2D Zhang-Suen Skeletonizing of Slice 90 




Figure 19: 2D Zhang-Suen Filling-Thinning-Contouring near Fovea 
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A Proof-of- Concept Code Excerpts 

A.l AScanLoader 

package marf . apps . oct . OCTMARF . Storage ; 

import Java. io .Buff eredReader ; 
import Java. io . File ; 
import Java. io . FileReader ; 
import Java. io . lOException; 

import marf . apps . oct . OCTMARF. framework. processors . ZhangSuenTransf orm; 
import marf .util . Arrays ; 

/** 

* Loads and parses A-scans from .txt files into internal data structures. 
* 

* Oauthor Serguei Mokhov 

* Osince 0.0.1, October 30, 2011 

* aversion $Id: AScanLoader. Java, v 1.8 2012/01/02 08:06:32 mokhov Exp $ 
*/ 

public class AScanLoader 
{ 

/** 

* 3D image volume storage. 
*/ 

protected Volume oVolume = new Volume (); 
/** 
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* Current slice to process. 

* Usually during loading. 
*/ 

protected int iScanSlice = 0; 
/** 

* Empty loader. 
*/ 

public AScanLoaderO 

{ 

} 

/** 

* Should the scan fail to load. 

* aparam pstrSampleFilename 

* Sthrows lOException 

* Ssee {aiink #load (String)} 
*/ 

public AScanLoader (String pstrSampleFilename) 

throws lOException 

{ 

load(pstrSampleFilename) ; 

} 

/** 

* Loads the current specified slice. 

* Sparam pstrSampleFilename 

* Sthrows lOException 
*/ 

public void load(String pstrSampleFilename) 

throws lOException 

{ 

Buf f eredReader oBFR = new Buff eredReader (new FileReader (pstrSampleFilename) ) ; 
String strDataLine = oBFR . readLine ( ) ; 

int iAScanLine = 0; 

while (strDataLine != null) 

//System . out .printlnC'Data line: [" + StrDataLine + "] ") ; 

String [] astrPixels = strDataLine . split (" "); 

//System. out. printlnC'Data line element count: [" + astrPixels . length + "]"); 

//this . astrSDScan [this . iScanSlice] [iAScanLine] = astrPixels; 
//Arrays. copy (this. as3DScan [this. iScanSlice] [iAScanLine] , astrPixels) ; 

//Arrays. copy (this. adSDScan [this. iScanSlice] [iAScanLine] , this. asSDScan [this. iScanSlice] [iAScanLine] ) ; 

this. oVolume.astrSDScan [this. iScanSlice] [iAScanLine] = astrPixels; 

Arrays. copy (this. oVolume.as3DScan [this. iScanSlice] [iAScanLine] , astrPixels) ; 

Arrays. copy (this. oVolume.adSDScan [this. iScanSlice] [iAScanLine] , this .oVolume . asSDScan [this . iScanSlice] [iAScanLine] ) ; 
iAScanLine++; 

StrDataLine = oBFR. readLine () ; 
//strDataLine = null; 

} 

oBFR.closeO ; 

// Move to the next slice for (potential) loading 
this . iScanSlice++ ; 

} 

/** 

* Sreturn 
*/ 

public short [][][] getShort3DScan() 
{ 
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return this . oVolime . asSDScan ; 

} 

/** 
* Oreturn 
*/ 

public double[][][] getDoubleSDScanO 
return this . oVolume . adSDScan ; 



/** 
* Sreturn 
*/ 

public String[][][] getString3DScan() 
return this . oVolume . astrSDScan; 



public Volume getVolumeO 
return this . oVolume; 



public int sizelnSlices () 
return this . iScanSlice + 1; 



// EOF 



A. 2 Volume 

package marf . apps . oct . OCTMARF . Storage ; 

import java.io.Serializable; 

/** 

* TODO: document. 
* 

* 

* $Id: Volume. Java, V 1.1 2011/12/18 07:17:34 mokhov Exp $ 
* 

* Sauthor serguei 

* Ssince 
*/ 

public class Volume 

implements Serializable, Cloneable 

{ 

/** 
* 

*/ 

private static final long serial VersionUID = -6995648761632634728L; 
/** 

* A-scan's line. 
*/ 

public static final int DEFAULT_SCAN_HEIGHT = 300; 
/** 

* Normal scan dimension. 
*/ 

public static final int DEFAULT_SCAN_WIDTH = 480; 
/** 

* Number of scan slices. 
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*/ 

public static final int DEFAULT_SCAN_DEPTH = 100; 

protected short [] [] [] asSDScan = new short [DEFAULT_SCAN_DEPTH] [DEFAULT_SCAN_W1DTH] [DEFAULT_SCAN_HEIGHT] ; 
protected double [] [] [] adSDScan = new double [DEFAULT_SCAN_DEPTH] [DEFAULT_SCAN_WIDTH] [DEFAULT_SCAN_HE1GHT] ; 
protected String [] [] [] astrSDScan = new String [DEFAULT.SCAN.DEPTH] [DEFAULT.SCAN.WIDTH] [DEFAULT.SCAN.HEIGHT] ; 



private boolean bEmpty = true; 

public VolumeO 
{ 

super ; 

} 

/** 

* Sreturn the asSDScan 
*/ 

public short [][][] getShortSDScanO 
{ 

return asSDScan; 

} 

/** 

* laparam asSdScan the asSDScan to set 
*/ 

public void setShort3DScan(short □ □ □ asSdScan) { 
asSDScan = asSdScan; 

} 

/** 

* Sreturn the adSDScan 
*/ 

public double [][] [] getDoubleSDScanO { 
return adSDScan; 

} 

/** 

* Oparam adSdScan the adSDScan to set 
*/ 

public void setDoubleSDScan (double [][] [] adSdScan) { 
adSDScan = adSdScan; 

} 

/** 

* Qreturn the astrSDScan 
*/ 

public String [] [] [] getStringSDScanO { 
return astrSDScan; 

} 

/** 

* Oparam astrSdScan the astrSDScan to set 
*/ 

public void setStringSDScan(String [] [] [] astrSdScan) { 
astrSDScan = astrSdScan; 

} 

public boolean isEmptyO 

{ 

return this . bEmpty ; 

} 

/* (non-Javadoc) 

* ®see java.lang.Object#clone() 
*/ 

QOverride 

public Object cloneO 
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//throws CloneNotSupportedException 
{ 

Volume oCopy = null; 

try 
{ 

oCopy = (Volume) super. clone () ; 

if (this.adSDScan != null) 
{ 

oCopy . adSDScan = this. adSDScan. clone ; 

} 

if (this.as3DScan != null) 
{ 

oCopy.as3DScan = this. asSDScan. clone ; 

} 

if (this.astrSDScan != null) 
{ 

oCopy . astrSDScan = this . astrSDScan. cloneO ; 

} 

} 

catchCCloneHotSupportedException e) 
{ 

System. err .printlnC "Failed to clone volume."); 
e.printStackTrace (System. err) ; 

} 

return oCopy; 

} 

/* (non-Javadoc) 

* ®see Java. lang. Qbject#f inalizeO 
*/ 

aOverride 

protected void finalizeO 

throws Throwable 

{ 

this . adSDScan = null; 
this . as3DScan = null; 
this.astrSDScan = null; 
super .finalizeO ; 

} 

/* (non-Javadoc) 

* ®see java.lang.Object#toString() 
*/ 

SOverride 

public String toStringO 
{ 

return "}: empty: " + this.bEmpty; 

} 



// EOF 



A.3 OCTMARFApp 

package marf . apps . oct . OCTMARF ; 

import Java. io .File ; 
import java.util.ArrayList; 

import marf . MARF ; 

import marf . Classificat ion. Distance .HammingDistance; 
import marf . Storage .ModuleParams ; 
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import marf . Storage .TrainingSet ; 

import marf . apps . oct . OCTMARF. Storage . AScanLoader ; 

import marf .apps . oct . OCTMARF. Storage . AnyToWAVLoader; 

import marf .apps . oct . OCTMARF. Storage . PPMDumper; 

import marf . apps . oct . OCTMARF . Storage . RawSampleLoader ; 

import marf . apps .oct . OCTMARF. framework. preprocessing. Binarizer; 

import marf . apps . oct . OCTMARF . framework . processors . ZhangSuenTransf orm; 

import marf .util. Debug; 

import marf . util . MARFException ; 

import marf .util .Not ImplementedExcept ion; 

import marf .util. OptionProcessor; 

/** 

* <p> 

* A an OCT image processing testbed platform for various algorithms. 

* </p> 
* 

* Sauthor Serguei Mokhov 
* 

* aversion $Id: OCTMARF App.java.v 1.6 2012/01/02 08:06:32 mokhov Exp $ 

* asince October 2011 
*/ 

public class OCTMARFApp 
{ 

/* 

* Apps. Versioning 

*/ 

/** 

* Current major version of the application. 
*/ 

public static final int MAJOR. VERSION = 0; 
/** 

* Current minor version of the application. 

*/ 

public static final int MINOR. VERSION = 0; 
/** 

* Current revision of the application. 
*/ 

public static final int REVISION = 1; 
/* 

* 

* Defaults 

* 

*/ 

/** 

* Default file name extension for sample files; case-insensitive. 
*/ 

public static final String DEFAULT.SAMPLE.FILE.EXTENSION = "txt"; 
/* 

* 

* Major and Hisc Options Enumeration 
« 

*/ 
/** 

* Numeric equivalent of the option <code> — train</code> . 
*/ 

public static final int OPT.TRAIN = 0; 
/** 

* Numeric equivalent of the option <code> — ident</code>. 
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*/ 

public static final int OPT.IDENT = 1; 
/** 

* Numeric equivalent of the option <code> — stats</code>. 
*/ 

public static final int OPT.STATS = 2; 
/** 

* Numeric equivalent of the option <code> — reset</code>. 
*/ 

public static final int OPT.RESET = 3; 
/** 

* Numeric equivalent of the option <code> — version</code> . 

*/ 

public static final int OPT_ VERSION = 4; 
/** 

* Numeric equivalent of the option <code> — help</code>. 
*/ 

public static final int OPT_HELP_LONG = 5; 
/** 

* Numeric equivalent of the option <code>-h</code> . 
*/ 

public static final int OPT_HELP_SHORT = 6; 
/** 

* Numeric equivalent of the option <code>-debug</code>. 
*/ 

public static final int OPT.DEBUG = 7; 
/** 

* Numeric equivalent of the option <code>-graph</code> . 
*/ 

public static final int OPT.GRAPH = 8; 

/** 

* Numeric equivalent of the option <code>-spectrogram</code> . 
*/ 

public static final int OPT.SPECTROGRAM = 9; 
/** 

* Numeric equivalent of the option <code>&lt ; writer ID></code>. 
*/ 

public static final int OPT_EXPECTED_ID = 10; 
/** 

* Numeric equivalent of the option <code> — batch-ident</code> . 
*/ 

public static final int OPT_BATCH_IDENT = 11; 
/** 

* Numeric equivalent of the option <code> — single-train</code> . 
*/ 

public static final int 0PT_SINGLE_TRA1N = 12; 
/** 

* Numeric equivalent of the option 

* <code>&lt ; sample-f ile-or-directory-name&gt ; </code> . 
*/ 

public static final int OPT_DIR_OR_FILE = 13; 
/** 

* Numeric equivalent of the option <code> — best-score</code>. 
*/ 

public static final int OPT_BEST_SCORE = 14; 
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/** 

* Numeric equivalent of the option <code> — gui</code> to have a GUI. 
*/ 

public static final int OPT.GUI = 15; 
/** 

* Numeric equivalent of the option <code>-noise</code> to enable noise 

* removal. 
*/ 

public static final int OPT_NOISE_REMOVAL = 16; 
/** 

* Numeric equivalent of the option <code>-silence</code> to enable silence 

* removal. 

*/ 

public static final int OPT_SILENCE_REMOVAL = 17; 
/** 

* Numeric equivalent of the option <code> — show-training-set</code> to 

* display text representation of a training set file. 
*/ 

public static final int OPT_SHOW_TRAINING_SET_DATA = 18; 

public static final int OPT_NO_PREPREP = 20; 

public static final int OPT.SILENCE.NOISE.REMOVAL = 27; 

public static final int OPT_TRAINING_MODEL_MEAN_CLUSTER = 29; 
public static final int OPT_TRAINING_MODEL_MEDIAN_CLUSTER = 30; 
public static final int OPT_TRAINING_MODEL_FEATURE_SET = 31; 

/* 

* Some defaults 
*/ 

public static final Integer DEFAULT_SEVERITY = 2; 
public static final Double DEFAULT_PROBABILITY = 0.0; 

public static final Double DISTANCE_THRESHOLD =0.1; 
public static final Double SIMILARITY_THRESHOLD = 0.001; 
public static final Double PROBABILITY.THRESHOLD =0.9; 

/* 

* 

* State Data Structures 
« 

*/ 
/** 

* Instance of the option processing utility. 
*/ 

protected static OptionProcessor soGetOpt = new OptionProcessor () ; 
/* 

* Static State Init 

*/ 

static 
{ 

// Main options 

soGetOpt.addValidOption(DPT_TRAIN, " — train") ; 
soGetDpt.addValidOption(OPT_SINGLE_TRAIN, " — single-train") ; 
soGetOpt.addValidOption(QPT_IDENT, " — Ident") ; 
soGetOpt.addValidOption(OPT_BATCH_IDENT, "— batch-ident") ; 
soGetOpt.addValidOption(OPT_STATS, "—stats") ; 
soGetOpt .addValidOption(OPT_BEST_SCORE, " — best-score") ; 

soGetOpt . addValidOpt ion (OPT_SHOW_TRAINING_SET_DATA , " —show-training-set " ) ; 



34 



soGetOpt.addValidOptioii(OPT_VERSION, " — version") ; 
soGetDpt.addValidOptioii(OPT_RESET, " — reset") ; 
soGetQpt.addValidOptioii(OPT_HELP_LONG, "—help") ; 
soGetOpt.addValidOptioii(OPT_HELP_SHORT, "-h") ; 
soGetOpt . addValidOptionCOPT.DEBUG, "—debug") ; 

// Sample Loading 

soGetOpt.addValidOptionCMARF.TEXT, "-text") ; 



// Preprocessing 

soGetOpt.addValidOption(OPT_HOISE_REMOVAL, "-noise"); 
soGetOpt.addValidOption(OPT_SILENCE_REMOVAL, "-silence") ; 
soGetOpt.addValidOption(OPT_SILENCE_KQISE_REMOVAL, "-silence-noise"); 



soGetDpt.addValidOption(MARF. DUMMY, "-norm") ; 

soGetOpt.addValidOption(MARF.HIGH_FREQUENCY_BOOST_FFT_FILTER, "-boost") ; 
soGetOpt.addValidOption(MARF.HIGH_PASS_FFT_FILTER, "-high") ; 
soGetOpt.addValidOption(MARF.LOW_PASS_FFT_FILTER, "-low") ; 
soGetOpt.addValidOption(MARF.BANDPASS_FFT_FILTER, "-band") ; 
soGetOpt.addValidOption(MARF.BAND_STOP_FFT_FILTER, "-bandstop") ; 
soGetDpt.addValidOption(MARF.HIGH_PASS_BOOST_FILTER, "-highpassboost") ; 
soGetOpt.addValidOptionCMARF.RAW, "-raw") ; 
soGetOpt.addValidOptionCMARF.ENDPOINT, "-endp") ; 
soGetOpt . addValidOption(MARF .HIGH_PASS_CFE_FILTER, "-highcf e" ) ; 
soGetOpt.addValidOption(MARF.LOW_PASS_CFE_FILTER, "-lowcfe") ; 
soGetOpt . addValidOpt ion (MARF . BAND_PASS_CFE_FILTER , "-bandcf e " ) ; 
soGetOpt. addValidOption(MARF.BAND_STOP_CFE_FILTER, "-bandstopcf e") ; 



// Wavelets 

soGetOpt. addValidOption(MARF.SEPARABLE_DWT_FILTER, "-sdwt") ; 
soGetOpt. addValidOption (MARF. DUAL_DTREE_DWT_FILTER, "-dualtree") ; 
soGetOpt. addValidOption(MARF.DYADIC_DWT_FILTER, "-dyadic") ; 



// Feature extraction 

soGetOpt. addValidOption(MARF.FFT, "-fft") ; 
soGetOpt. addValidOption(MARF.LPC, "-Ipc") ; 

soGet Opt . addVal idOpt i on (MARF . RANDOM_FEATURE_EXTRACTI ON , " -r andf e " ) ; 

soGetOpt. addValidOption (MARF. MIN_MAX_AMPLITUDES, "-minmax") ; 

soGet Opt . addVal idOpt i on (MARF . FEATURE_EXTRACTION_ AGGREGATOR , " -aggr " ) ; 

soGetOpt. addValidOption(MARF.FO, "-fO") ; 

soGetOpt. addValidOption(MARF. SEGMENTATION, "-segm") ; 

soGetOpt. addValidOption(MARF.CEPSTRAL, "-cepstral") ; 



// Classification 

soGetOpt. addValidOption(MARF.NEURAL_NETWORK, "-nn") ; 
soGetOpt. addValidOption(MARF.EUCLIDEAN_DISTANCE, "-eucl") ; 
soGetOpt. addValidOption(MARF.CHEBYSHEV_DISTANCE, "-cheb") ; 
soGetOpt. addValidQption(MARF.MINKOWSKI_DISTANCE, "-mink") ; 
soGetOpt. addValidOption(MARF.MAHALANOBIS_DISTANCE, "-mah") ; 
soGetOpt. addValidOption(MARF.RANDOM_CLASSIFICATION, "-randcl") ; 
soGetOpt . addValidOption (MARF . DIFF.DISTANCE , " -dif f " ) ; 
soGetOpt. addValidOption(MARF.ZIPFS_LAW, "-zipf ") ; 
soGetOpt. addValidOption(MARF. MARKOV, "-markov") ; 
soGetOpt. addValidOption(MARF.HAMMING_DISTANCE, "-hamming") ; 
soGetOpt . addValidOption (MARF . COSINE.SIMILARITY.MEASURE , " -cos " ) ; 

// Graphical 

soGetOpt. addValidOption(OPT_SPECTROGRAM, "-spectrogram") ; 
soGetOpt. addValidOption(OPT_GRAPH, "-graph") ; 
soGetOpt. addValidOption(OPT_GUI, " — gui") ; 

// NLP 

soGetOpt .addValidOption (MARF. EStatisticalEstimators.MLE, "-mle") ; 

soGetOpt .addValidOption (MARF. EStatisticalEstimators.ADD.OWE, "-add-one") ; 

soGetOpt .addValidOption (MARF. EStatisticalEstimators.ADD.DELTA, "-add-delta") ; 

soGetOpt .addValidOption (MARF. EStatisticalEstimators.ELE, "-ele") ; 

soGetOpt . addValidOption (MARF . EStatisticalEstimators . WITTEN.BELL , "-witten-bell" ) ; 
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soGetOpt .addValidOptioii(MARF.EStatisticalEstlmators.GOOD_TURING, "-good-turing") ; 

soGetOpt .addValidOptionCMARF.EStatisticalEstlmators.SLI, "-sll") ; 

soGetOpt .addValidOptionCMARF.EStatisticalEstimators.GLI, "-gll") ; 

soGetQpt .addValidOption(MARF.EStatisticalEstimators.KATZ_BACKQFF, "-backoff") ; 

soGetOpt. addValidOptionCMARF.WLP. STEMMING, "-stemming") ; 

soGetOpt. addValidOption(MARF.NLP.CHARACTER_MODE, "-char") ; 
soGetOpt . addValidOpt ion (MARF . NLP . WORD.MODE , " -word" ) ; 

soGetOpt. addValidQption(MARF.NLP.CASE_SENSITIVE, "-case") ; 
soGetOpt . addValidOption (MARF .NLP . PARSE.NUMBERS , "-num" ) ; 
soGetOpt . addValidOption (MARF . NLP . PARSE.QUOTED.STRINGS , " -quote" ) ; 
soGetOpt . addValidOption (MARF . NLP . PARSE_ENDS_DF_SENTENCE , " -eos " ) ; 
soGetOpt . addValidOption (MARF . NLP . RAW_ZIPFS_LAW_DUMP , " -nolog" ) ; 

soGetOpt . addValidOption (MARF . ENgramModels . UNIGRAM , "-unigram" ) ; 
soGetOpt. addValidOption(MARF.ENgramModels.BIGRAM, "-bigram") ; 
soGetOpt . addValidOption (MARF . ENgramModels . TRIGRAM , "-trigram" ) ; 

soGetOpt. addValidOption(MARF. NLP. TRAIN, " — train-nip") ; 
soGetOpt. addValidOption(MARF. NLP. CLASSIFY, "— ident-nlp") ; 

soGetOpt. addValidOption(MARF. NLP. ZIPFS.LAW.CHEAT, "—cheat") ; 

soGetOpt. addValidOption(MARF. NLP. INTERACTIVE, "-interactive") ; 

soGetOpt . addValidOption(OPT_NO_PREPREP, "-nopreprep") ; 



public static class Configuration 
{ 

String strFilePrefix = ""; 
boolean bBinarize = false; 
boolean bZSSkeletonlzing = false; 
boolean bZSFilling = false; 
boolean bZSThinning = false; 
boolean bZSCountouring = false; 

int iRepeat = 1 ; 

public Conf igurationO 
{ 



public String toStringO 
{ 

return strFilePrefix + " : " + bBinarize + " : " + bZSSkeletonlzing + " : " + bZSFilling + " : " + bZSThinning + " : " + bZSCounto 

} 

} 



/** 

* Main body. 
* 

* Sparam argv command-line arguments 
*/ 

public static final void main(String[] argv) 
{ 

try 
{ 

// Since some new API is always introduced. . . 
validateVersionsO ; 

// Parse extra arguments 

int iValidOptions = soGetOpt .parse(argv) ; 



// Turn on command-line supplied debugging 
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Debug. enableDebug(soGetOpt .isActiveOptioii(DPT_DEBUG)) ; 
Debug. debug ("Opt ion set parsed: " + soGetOpt) ; 

if (iValidOptions == 0) 
{ 

throw new ExceptionC'No valid options found: " + soGetOpt) ; 

} 

switch ( soGetOpt . get InvalidOpt ions O.sizeO) 
{ 

case 0: 
{ 

break; 

} 

/* 

* In the case of a single invalid option always assume it 

* is either a filename or a directory name for the — ident 

* or — train options. 
*/ 

case 1: 
{ 

soGetOpt. addActiveOption(OPT_DlR_QR_FILE, soGetOpt .getlnvalidOptions () .f irstElementO .toStringO) ; 

soGetOpt .getlnvalldOptionsO . clear () ; 

break; 

} 

default : 
{ 

throw new ExceptionCUnrecognized options found: " + soGetOpt .getlnvalidOptionsO ) ; 

} 



setDef aultConf igO ; 
setCustomConf igO ; 

// Set misc configuration 

MARF . setDumpSpectrogramCsoGetOpt . isActiveOption(OPT_SPECTROGRAM) ) ; 
MARF. setDumpWaveGraphCsoGetOpt . isActiveOption(OPT_GRAPH) ) ; 

Debug. debugC'Option set: " + soGetOpt); 

if (soGetOpt . isActiveOptionCMARF . TEXT) ) 
{ 

MARF . setSampleFormat (MARF . TEXT) ; 

MARF . setSampleLoaderPluginClass (AnyToWAVLoader . class . getName ( ) ) ; 
System. out. printlnC'Using text loader"); 

} 

else if (soGetOpt. isActiveOption (MARF. TEXT +2)) 
{ 

MARF . setSampleFormat (MARF . CUSTOM) ; 
// XXX: 

MARF . SetSampleLoaderPluginClass (RawSampleLoader . class . getName ) ; 



int iMainOption = soGetOpt .getOption(argv[0] ) ; 

switch ( iMainOpt ion) 
{ 

/* 

* 

* Identification 
* 

*/ 

// Single case 

case OPT.IDENT: 

case OPT_BATCH_IDEHT: 
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case MARF.NLP. CLASSIFY: 
{ 

ArrayList<Conf iguration> oConfigs = new ArrayList<Conf iguration>() ; 

Configuration oConf ig = new Conf igurationO ; 

oConf ig.bBinarize = false; 

oConf ig.bZSCoimtouring = false; 

oConf ig.bZSFilling = false; 

oConf ig.bZSThimiing = false; 

oConf ig.bZSSkeletonizing = true; 

oConf ig. strFilePref ix = "zsuen-skel" ; 

oConf ig. iRepeat = 1; 

oConfigs .add(oConfig) ; 

oConf ig = new Conf igurationO ; 
oConf ig.bBinarize = true; 
oConf ig.bZSCountouring = false; 
oConf ig.bZSFilling = false; 
oConf ig.bZSThinning = false; 
oConf ig.bZSSkeletonizing = false; 
oConfig. StrFilePref ix = "simple-bin"; 
oConfig. IRepeat = 1; 
oConfigs . add(oConfig) ; 

oConfig = new Conf igurationO ; 
oConf ig.bBinarize = false; 
oConf ig.bZSCountouring = true; 
oConf ig.bZSFilling = true; 
oConf ig.bZSThinning = true; 
oConf ig.bZSSkeletonizing = false; 
oConfig. StrFilePref ix = "zsuen-ftc"; 
oConfig. iRepeat = 1; 
oConfigs .add(oConfig) ; 

oConfig = new Conf igurationO ; 
oConf ig.bBinarize = false; 
oConf ig.bZSCountouring = false; 
oConf ig.bZSFilling = true; 
oConf ig.bZSThinning = true; 
oConf ig.bZSSkeletonizing = true; 
oConfig. StrFilePref ix = "zsuen-f ts-3" ; 
oConf ig . iRepeat = 3 ; 
oConfigs .add(oConfig) ; 

/* 

oConfig = new Conf igurationO ; 

oConf ig.bBinarize = false; 

oConf ig.bZSCountouring = false; 

oConf ig.bZSFilling = true; 

oConf ig.bZSThinning = false; 

oConf ig.bZSSkeletonizing = true; 

oConfig. StrFilePref ix = "zsuen-f s-3" ; 

oConfig. iRepeat = 3; 

oConfigs .add(oConfig) ; 

*/ 

for (Conf iguration oRunConfig: oConfigs) 
{ 

ident (oRunConf ig) ; 

} 



break; 

} 



/* 

* 

* Training 
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*/ 

// Add a single sample to the training set 

case 0PT_S1NGLE_TRAIN : 

{ 

throw new NotlmplementedExceptionO ; 

} 

// Train on a directory of files 
case OPT.TRAIK: 
case MARF.KLP.TRAIH: 
{ 

throw new NotlmplementedExceptionO; 

} 

/* 

* 

* Training Set Data Display 
* 

*/ 

case OPT_SHOW_TRAINING_SET_DATA: 
{ 

System. out .println 
( 

marf .Classification. Classification. loadTrainingSet 
( 

marf . Storage . IStorageManager . DUMP_GZIP_BINARY , 
soGetOpt.getOption(OPT_DIR_OR_FILE) 

) 

); 

break; 

} 

/* 
* 

* Stats 
* 

*/ 

case OPT_STATS: 
{ 

throw new NotlmplementedExceptionO; 

} 

/* 

* Best Result with Stats 
*/ 

case OPT_BEST_SCORE: 
{ 

throw new NotlmplementedExceptionO; 

} 
/* 

* Reset Stats 
*/ 

case OPT.RESET: 
{ 

throw new NotlmplementedExceptionO; 

} 

/* 

* Versioning 
*/ 

case OPT.VERSION: 
{ 

System. out. printlnC'OCTMARF Application, v." + getVersionO ) ; 
System. out. printlnC'Using MARF, v." + MARF. get Ver s ionO ) ; 
validateVersionsO ; 
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break; 

} 

/* 

* Help 
*/ 

case OPT.HELP.SHORT: 
case OPT_HELP_LONG: 
{ 

usage () ; 
break; 

} 

/* 

* Invalid major option 
*/ 

default : 
{ 

throw new ExceptionC'Unrecognized option: " + argv[0]); 

} 

y // major option switch 
} // try 

/* 

* No arguments have been specified 
*/ 

catch(ArrayIndexOutOf BoundsException e) 
{ 

System. err. printlnC'Ko arguments have been specified."); 
System. err .printlnCe .getMessageO ) ; 
e.printStackTrace (System. err) ; 

usage () ; 

} 

/* 

* MARF-specif ic errors 
*/ 

catch(MARFException e) 
{ 

System, err .printlnCe .getMessageO ) ; 
e.printStackTrace (System. err) ; 

} 

/* 

* Invalid option and/or option argument 
*/ 

catch(Exception e) 
{ 

System. err .println(e. getMessageO ) ; 
e . printStackTrace (System . err) ; 
usage O ; 

} 

/* 

* Regardless whatever happens, close the db connection. 
*/ 

finally 
{ 

try 
{ 

Debug. debugC'Closing DB connection..."); 

} 

catch (Except ion e) 
{ 

Debug. debug("Closing DB connection failed: " + e.getMessageO) ; 
e.printStackTrace (System. err) ; 
System. exit(-l) ; 
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} 

} 

} 

/** 

* Computes duration and outputs it as a string 

* 8param plStartTime 

* Sparam plEndTime 

* Sreturn String-formatted duration 
*/ 

private static String getDurat ion (long plStartTime, long plEndTime) 
{ 

// XXX: constify 

int iOneSecInMs = 1000; 

int iOneMinlnMs = 60 * iOneSecInMs; 

int lOneHourlnMs = 60 * iOneMinlnMs; 
int iOneDaylnMs = 24 * iOneHourlnMs ; 

long IDuration = plEndTime - plStartTime; 

long IDays = IDuration / iOneDaylnMs; 

long IModTimeLeft = IDuration '/, iOneDaylnMs; 

long IHours = IModTimeLeft / iOneHourlnMs; 
IModTimeLeft '/,= iOneHourlnMs; 

long IMinutes = IModTimeLeft / iOneMinlnMs; 
IModTimeLeft '/.= iOneMinlnMs; 

long ISeconds = IModTimeLeft / IOneSecInMs; 
IModTimeLeft '/.= iOneSecInMs; 

long IMilliSeconds = IModTimeLeft; 

StringBuffer oDuration = new StringBuf f er() 

. appendClDays + "d:" + IHours + "h:" + IMinutes + "m:" + ISeconds +"s:" + IMilliSeconds +"ms:" + IDuration + "ms") 



return oDuration. toStringO ; 

} 

/** 

* Generic ident routine. 

* ®throws Exception 
*/ 

public static final void ident (Configuration poConfig) 

throws Exception 

{ 

long ITimeStart = System. currentTimeMillis ; 
{ 

AScanLoader oLoader = new AScanLoader () ; 
// Process all .txt files in the directory 

File[] aoSampleFiles = new File("resources/data") .listFilesO ; 

Binarizer oBinarizer = new BinarizerO ; 

ZhangSuenTransf orm oTransf orm = new ZhangSuenTransf orm() ; 

for(int i = 0; i < aoSampleFiles . length; i++) 
{ 

String strFilename = aoSampleFiles [i] .getPathO ; 

if (aoSampleFiles [i] .IsFlleO && strFilename. toLowerCaseO . endsWith(DEFAULT_SAMPLE_FILE_EXTENS10N) ) 
{ 

long ISliceTimeStart = System. currentTimeMillis ; 
// Load them one by one 

System. out .printlnC'Loading " + strFilename + " into slice " + oLoader.sizelnSlicesO + " of 3D volume."); 
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oLoader.load(strFilename) ; 

// Trial processing 
{ 

PPMDumper oDumper = new PPMDumper(poConf ig.strFilePrefix + "." + i) ; 

if (poConf ig.bBinarize == true) 
{ 

oBinarizer.preprocess2D(oLoader. get Volume , i) ; 

} 

if 
( 

poConf ig.bZSCountouring == true 
I I 

poConf ig.bZSFilling == true 
I I 

poConf ig.bZSSkeletonizing == true 
I I 

poConf ig.bZSThinning == true 

) 
{ 

oTransf orm . setCountour ing (poConf ig . bZSCountour ing) ; 
oTransform. setFillingCpoConf ig.bZSFilling) ; 
oTransf orm. setSkeletonizing (poConf ig.bZSSkeletonizing) ; 
oTransf orm . setThinning (poConf ig . bZSThinning) ; 

for(int r = 0; r < poConf ig. iRepeat ; r++) 
{ 

oTransform. detect2D(oLoader.getVolume() , i) ; 

} 

} 

System. out. printC'Dumping " + oDumper . getFilename () + "... "); 

for(double[] adScanLine: oLoader .getDouble3DScan() [i] ) 

{ 

oDumper .addToData(adScanLine) ; 

} 

// Save as a PPM 
oDumper . dump ( ) ; 
System. out. printlnC'Done. ") ; 
} // end trial processing 

long ISliceTimeFinish = System. currentTimeMillis () ; 

System . out . println 

( 

"Slice processing and saving time: " + 
getDurationdSliceTimeStart , ISliceTimeFinish) 

); 

} // if 

} 

} 

long ITimeFinish = System. currentTimeMillis () ; 

System. out .println 

( 

"Total time volume took: " + getDuration(lTimeStart , ITimeFinish) 
+ " for configuration: " + poConfig 

); 

} 

/** 

* Displays application's usage information and exits. 
*/ 

private static final void usageO 
{ 

System . out . println 
( 



'Usage : \n" 

Java OCTMARFApp 



— train <samples-dir> [options] 
— single-train <sample> [options] 
— ident <sample> [options] 
— batch-ident <sainples-dir> [options 
—train-nip [ —debug ] [ OPTIONS ] 
— ident-nlp [ —debug ] [ OPTIONS ] 
— debug 
— gui 

— stats= [per-conf ig I per-wr iter I both] 
— best-score 
— reset 
— version 
—help I -h 



- train mode\n" 

- add a single sample to the training set\r 

- identification modeXn" 
] — batch identification mode\n" 
<language> <corpus-f ile>\n" 
foo <bar|corpus-file>\n\n" 

— include verbose debug output\n" 

— use GUI as a user interfaceXn" 

— display stats (default is per-conf ig)\n" 

— display best classification result \n" 

— reset statsNn" 

— display version info\n" 

— display this help and exit\n\n" + 



Options (one or more of the following) : \n\n" + 



'Loaders : \n\n" 

H " -wav - assume WAVE files loading (default) \n" 

H " -text - assume loading of text samplesXn" 

I- " -tiff - assume loading of TIFF samples\n" 

V "\n" + 



'Preprocessing: \n\n" 



" -silence 


- remove silence (can be combined with any of the below) \n" 


" -noise 


- remove noise (can be combined with any of the below) \n" 


" -raw 


- no preprocessing\n" 


" —norm 


- use just normalization, no filtering\n" 


" -low 


- use low-pass EFT filter\n" 


" -high 


- use high-pass FFT filter\n" 


" -boost 


- use high-frequency-boost FFT preprocessor\n" 


" -band 


- use band-pass FFT filter\n" 


" -bandstop 


- use band-stop FFT filter\n" 


" -endp 


- use endpointing\n" 


" -lowcfe 


- use low-pass CFE filter\n" 


" -highcfe 


- use high-pass CFE filter\n" 


" -bandcfe 


- use band-pass CFE filter\n" 


" -bandstopcfe 


- use band-stop CFE filter \n" 


"\n" + 





'Feature Extraction: \n\n" 



-Ipc - use LPC\n" 

-fft - use FFT\n" 

-minmax - use Min/Max Amplitudes\n" 

-randfe - use random feature extraction\n" 

-aggr - use aggregated FFT+LPC feature extraction\n" 

-fO - use FO (pitch, or fundamental frequency; NOT IMPLEMENTED\n" 

-segm - use Segmentation (NOT IMPLEMENTED) \n" 

-cepstral - use Cepstral analysis (NOT IMPLEMENTED) \n" 
\n" + 



Classification: \n\n" 



" -ml 


- use 


Neural Network\n" 


" -cheb 


- use 


Chebyshev Distance\n" 


" -eucl 


- use 


Euclidean Distance\n" 


" -mink 


- use 


Minkowski Distance\n" 


" -diff 


- use 


Dif f-Distance\n" 


" -zipf 


- use 


Zipf's Law-based classif ier\n" 


" -randcl 


- use 


random classif ication\n" 


" -markov 


- use 


Hidden Markov Models (NOT IMPLEMENTED) \n" 


" -hamming 


- use 


Hamming Distance\n" 


" -cos 


- use 


Cosine Similarity Measure\n" 


"\n" + 







'NLP/Ngrams/Smoothing : \n\n" 

y " -interactive - interactive mode for classification instead of reading from a file\n" 
y " -char - use characters as n-grams (should always be present for this app)\n\n" 

V " -unigram - use UNIGRAM model\n" 
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+ " 


-bigram 


- use 


BIGRAM model\n" 


+ " 


-trigram 


- use 


TRIGRAM model\n\n" 


+ " 


-mle 


- use 


MLE\n" 


+ " 


-add-one 


- use 


Add-One smoothingXn" 


+ " 


-add-delta 


- use 


Add-Delta (ELE, d=0.5) smoothing\n 


+ " 


-witten-bell 


- use 


Witten-Bell smoothingXn" 


+ " 


-good-turing 


- use 


Good-Turing smoothingXn" + 



"Xn" 

+ " -spectrogram 
+ " -graph 
+ " <integer> 
+ "Xn" 



dump spectrogram image after feature extractionXn" 

dump wave graph before preprocessing and after feature extractionXn" 

expected subject IDXn" 



); 



System. exit (0) ; 



/** 

* Retrieves String representation of the application's version. 
* 

* ©return version String 
*/ 

public static final String getVersionO 
{ 

return MAJOR. VERSION + "." + MINOR. VERSION + "." + REVISION; 

} 

/** 

* Retrieves integer representation of the application's version. 
* 

* fireturn integer version 
*/ 

public static final int getlntVersionO 
{ 

return MAJOR. VERSION * 100 + MINOR. VERSION * 10 + REVISION; 

} 



/** 

* Makes sure the applications isn't run against older MARF version. Exits 

* with 1 if the MARF version is too old. 
*/ 

public static final void validateVersionsO 
{ 

if (MARF.getDoubleVersionO < (0 * 100 + 3 * 10 + + .6)) 
{ 

System . err . println 



( 



"Your MARF version (" + MARF.getVersionO 

+ ") is too old. This application requires 0.3.0.6 or above." 



System. exit(l) ; 



/** 

* Composes the current configuration of in a string form. 
* 

* ®param pstrArgv 

* set of configuration options passed through the command line; 

* can be null or empty. If latter is the case, MARF itself is 

* queried for its numerical set up inside. 
* 

* Sreturn the current configuration setup 
*/ 

public static final String getConf igStringCString [] pstrArgv) 
{ 
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// Store config and error/ successes for that config 

String strConfig = ""; 

if (pstrArgv ! = null kk pstrArgv . length > 2) 
{ 

// Get config from the command line 
for(int i = 2; i < pstrArgv. length; i++) 
{ 

// We do NOT want numerical IDs 

try 

{ 

Integer. parseint (pstrArgv [i]) ; 

} 

catch(NumberFormatException e) 
{ 

// Skip .xml file names 

if (pstrArgv [i] .matches (" . *\\ . xml") ) 

{ 

continue ; 

} 

StrConfig += pstrArgv [i] + " "; 

} 

} 

} 

else 
{ 

// Query MARF for it's current config 
StrConfig = MARF . getConf ig ( ) ; 

} 

return strConfig; 

} 

/** 

* Sets default MARF configuration parameters as normalization for 

* preprocessing, FFT for feature extraction, Euclidean distance for 

* training and classification with no spectrogram dumps and no debug 

* information, assuming WAVE file format. 
* 

* athrows MARFException 
*/ 

public static final void setDefaultConf ig() throws MARFException 
{ 

/* 

* Default MARF setup 
*/ 

MARF . setPreprocessingMethod(MARF . DUMMY) ; 

MARF. setFeatureExtractionMethod (MARF. FFT) ; 

MARF. setClassificationMethod (MARF. EUCLIDEAK.DISTANCE) ; 

MARF. setDumpSpectrogram(f alse) ; 

// MARF . setSampleFormat (MARF . WAV) ; 

// MARF . setSampleFormat (MARF . TEXT) ; 

MARF . setSampleFormat (MARF . CUSTOM) ; 

MARF . setSampleLoaderPluglnClass (AnyToWAVLoader. class .getName () ) ; 
//Debug. enableDebug(f alse) ; 

} 

/** 

* Customizes MARF's configuration based on the options. 
* 

* athrows MARFException 

* if some options are out of range 
*/ 

public static final void setCustomConf ig() 

throws MARFException 

{ 
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ModulePar ams oPar ams = new ModuleParams ( ) ; 

for 
( 

int iPreprocessingMethod = MARF.MIN_PREPROCESSING_METHDD; 
iPreprocessingMethod <= MARF.MAX.PREPROCESSING.METHOD; 
iPreprocessingMethod++ 

) 
{ 

if (soGetOpt . isActiveOption(iPreprocessingMethod) ) 
{ 

MARF.setPreprocessingMethod (iPreprocessingMethod) ; 

switch(iPreprocessingMethod) 
{ 

// Endpointing did not respond well to silence removal 
// and the definition of RAW assumes none of that is done. 
// Endpointing also has some extra params of its own. 
case MARF.RAW: 
case MARF.ENDPOIHT: 

case MARF. DUMMY: 

case MARF.HIGH_FREqUENCY_BOOST_FFT_FlLTER: 

case MARF.HIGH_PASS_FFT_FILTER: 

case MARF.L0W_PASS_FFT_F1LTER: 

case MARF.BANDPASS.FFT.FILTER: 

case MARF.BAND_STOP_FFT_FILTER: 

case MARF.HIGH_PASS_BOOST_FILTER: 

case MARF.BAND_PASS_CFE_FILTER: 

case MARF.BAND_STOP_CFE_FILTER: 

case MARF.LOW_PASS_CFE_FILTER: 

case MARF.HIGH_PASS_CFE_FILTER: 

{ 

// Kludge 

if (soGetOpt . isActiveOption(OPT_SILENCE_NOISE_REMOVAL) ) 
{ 

oParams . addPreprocessingParajii(true) ; 
oParams . addPreprocessingParam(true) ; 

} 

else 
{ 

// Normalization and filters seem to respond better 
// to silence removal. 

// The setting of the third protocol parameter 

// (silence threshold) 

//is yet to be implemented here. 

oParams . addPreprocessingParam(soGetOpt . isActiveOption(OPT_NOISE_REMOVAL) ) ; 
oParams . addPreprocessingParam(soGetOpt . isActiveOption(OPT_SILENCE_REMOVAL) ) ; 

} 

break; 

} 

default: 
{ 

assert false : "Not implemented valid preprocessing configuration parameter: " + iPreprocessingMethod; 

} 

} // switch 
break; 

} 

} 

for 
( 

int iFeatureExtractionMethod = MARF . MIN.FEATUREEXTRACTION.METHOD ; 
iFeatureExtractionMethod <= MARF.MAX.FEATUREEXTRACTION.METHOD; 
iFeatur eExtr act ionMethod++ 



46 



) 
{ 

if (soGetOpt . isActiveOption(iFeatureExtractionMethod) ) 
{ 

MARF . setFeatureExtractionMethod(iFeatureExtractionMethod) ; 

swit ch ( iFeatureExtr act ionMethod) 
{ 

case MARF.FFT: 
case MARF.LPC: 

case MARF . RANDOM.FEATURE.EXTRACTION : 

case MARF.MIN_MAX_AMPLITUDES: 

case MARF.FO: 

case MARF.CEPSTRAL: 

case MARF. SEGMENTATION: 

// For now do nothing; customize when these methods 

// become parameterizable . 

break; 

case MARF . FEATURE.EXTRACTIDN.AGGREGATOR: 
{ 

// For now aggregate FFT followed by LPC until 
// it becomes customizable 

oPairams.addFeatureExtractionParaiiiCnew Int eger (MARF. FFT) ) ; 

oParams . addFeatureExtractionParam(null) ; 

oParams . addFeatureExtractionParam(new Integer (MARF . LPC) ) ; 

oParams . addFeatureExtractionParam(null) ; 

break; 

} 

default: 

assert false; 
} // switch 

break; 

} 

} 

for 

( 

int iClassif icationMethod = MARF.MIN_CLASSIFICATION_METHOD; 
iClassificationMethod <= MARF.MAX.CLASSIFICATION.METHOD; 
iClassif icationMethod++ 

) 
{ 

if (soGetOpt. isActiveOption(iClassif icationMethod)) 
{ 

MARF.setClassificationMethodCiClassif icationMethod) ; 

swit ch(iClassif icationMethod) 
{ 

case MARF.NEURAL.NETWORK: 
{ 

// Dump/Restore Format of the TrainingSet 

oParams . addClassif icationParam(TrainingSet .DUMP_GZIP_BINARY) ; 

// Training Constant 

oParams. addClassif icationParam(0. 5) ; 

// Epoch number 

oParams . addClassif icationParam(20) ; 
// Min. error 

oParams . addClassif icationParam(0. 1) ; 



break; 

} 
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case MARF.HAMMING.DISTANCE: 

{ 

// Dump/Restore Format 

oParams . addClassif icatlonParam(TrainingSet .DUMP_GZIP_BINARY) ; 
// Type of hamming comparison 

oParams . addClassif icationParam(HammingDistance . STRICT_DOUBLE) ; 
break; 

} 

case MARF.MINKOWSKI.DISTANCE: 
{ 

// Dump/Restore Format 

oParams . addClassif icationParam(TrainingSet .DUMP_GZIP_BINARY) ; 

// Minkowski Factor 

oParams. addClassif icationParam(6.0) ; 

break; 

} 



case 


MARF. 


, EUCLIDEAN_DISTANCE : 


case 


MARF. 


, CHEBYSHEV_DISTANCE : 


case 


MARF. 


, MAHALANOBIS.DISTANCE : 


case 


MARF. 


, RANDOM.CLASSIFICATION : 


case 


MARF. 


.DIFF.DISTANCE: 


case 


MARF. 


.MARKOV: 


case 


MARF. 


.ZIPFS.LAW: 


case 


MARF. 


. COSINE_SIMILARITY_MEASURE 


{ 

// 


Dump/Restore Format 



oParams . addClassif icationParam(TrainingSet .DUMP_GZIP_BINARY) ; 

// For now do nothing; customize when these methods 

// become parameterizable . 

break; 

} 

default : 

assert false : "Unrecognized classification module"; 

} // switch 

// Method is found, break out of the look up loop 
break; 

} 

} 

// Assign meaningful params only 

if (oParams. size () > 0) 

{ 

MARF.setModuleParams(oParams) ; 

} 

} 

} 

// EOF 



A.4 PPMDumper 

package marf . apps . oct . OCTMARF . Storage ; 

import Java. io .Buff eredOutputStream; 
import java.io.DataOutputStream; 
import Java. io .FileOutputStream; 
import java.util. Vector; 
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import marf . Storage . StorageException; 
import marf . util . Debug ; 

/** 

* <p>Dumps an image to a PPM file.</p> 
* 

* 9author Serguei Mokliov 

* aversion $Id: PPMDumper . java.v 1.5 2012/01/02 08:06:33 mokhov Exp $ 

* Ssince 0.0.1 
*/ 

public class PPMDumper 
{ 

/** 

* The data vector. 
*/ 

protected Vector<Double [] > oData = null; 
/** 

* Current minimum. 
*/ 

protected double dMin = 0.0; 
/** 

* Current maximum. 
*/ 

protected double dMax = 0.0; 
/** 

* To differentiate file names based on the feature extraction method name. 
*/ 

protected String strMethod = ""; 
/** 

* Constructor. 
*/ 

public PPMDumperO 
{ 

this.oData = new Vector<Double □ >() ; 

} 

/** 

* Constructor with a feature extraction method name. 

* (Sparam pstrMethodWame String representing FE module name 
*/ 

public PPMDumper (String pstrMethodName) 
{ 

thisO ; 

this . StrMethod = pstrMethodName; 

} 

/** 

* Dumps image. 

* ©throws StorageException 
*/ 

public final void dumpO 
throws StorageException 
{ 

try 
{ 

Debug. debug ("Dumping image " + this . strMethod + ".ppm"); 

Debug. debugC .dump () - data size in vectors: " + this.oData.sizeO) ; 

FileOutputStream oFOS = null; 
DataOutputStream oOutFile = null; 

// oFOS = new File0utputStream("58. " + this . strMethod + ".ppm"); 

//oFOS = new File0utputStream("58.ppm") ; 
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oFOS = new FileOutputStreamCthis . strMethod + ".ppm"); 
//oOutFile = new DataOutputStream(oFOS) ; 

oOutFile = new DataOutputStream(new BufferedOutputStream(oFOS)) ; 
// Output PPM header 

/* 

man ppm: 

- A "magic number" for identifying the file type. A pgm file's magic number is the two characters "P6". 

- Whitespace (blanks, TABs, CRs, LFs) . 

- A width, formatted as ASCII characters in decimal. 

- Whitespace. 

- A height, again in ASCII decimal. 

- Whitespace. 

- The maximum color value (Maxval) , again in ASCII decimal. Must be less than 65536. 

- Newline or other single whitespace character. 

- A raster of Width * Height pixels, proceeding through the image in normal English reading order. Each pixel is a 
green, and blue samples, in that order. Each sample is represented in pure binary by either 1 or 2 bytes. If 
than 256, it is 1 byte. Otherwise, it is 2 bytes. The most significant byte is first. 

*/ 

// Output data 

// Make max be at 757, 
this.dMax = 255; 

int iWidth = this. oData. size () ; 

int iHeight = this. oData. element At (0) .length; 

oOutFile . writeBytes 
( 

"P6\n" 

+ iWidth + "\n" 
+ iHeight + "\n" 
+ "255\n" 

); 

for(int i = iHeight - 1; i >= 0; i — ) 
{ 

forCint j = 0; j < iWidth; 

Double □ adData = this.oData.elementAt(j) ; 

// Colors, RGB, hence 3 
for(int m = 0; m < 3; m++) 

{ 

oOutFile. writeByte (adData [i] .byteValueO) ; 

} 

} 

} 

oOutFile. flushO ; 

Debug . debug 
( 

"Done dumping image " + 
this . StrMethod + 

".ppm [" + (this . oData. sizeO * iHeight * iWidth) + " bytes]" 

); 

} 

catch(Exception e) 
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{ 

throw new StorageException(e) ; 

} 

} 

public void addToData (double [] padData) 
{ 

Double [] oNewData = new Double [padData. length] ; 
//Arrays. copyCoNewData, paData) ; 

for(int i = 0; i < padData. length; i++) 
{ 

oNewData[i] = padData [i] ; 

} 

this . oData. add(oNewData) ; 

} 

public void addToData(Double □ paData) 
{ 

this . oData. add (paData) ; 

} 

/** 

* Oreturn the oData 
*/ 

public Vector<Double □ > getDataO 
{ 

return this.oData; 

} 

/** 

* Sparam oData the oData to set 
*/ 

public void setData(Vector<Double[]> poData) 
{ 

this.oData = poData; 

} 

public String getFilename () 
{ 

return this . strMethod ; 

} 



// EOF 
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